// EnumerateTIsoforms.cpp : Defines the entry point for the console application.
//
#include "IsoformsEnumerator.h"
#include "IsoformsSelecter.h"
#include "outputIsoforms.h"
#include "Parsers.h"
#include "Junctions.h"
#include "chromosomeNTdata.h"
#include "stdafx.h"
#include <iostream>
#include <vector>
using namespace std;
const int MAX_EXONS_NO = 30;
const int MAX_NO_ISOFORM_CANDIDATE = 20;
void print_info(void)
{
    printf("The only arguement should be the <gene name>.\n");
    printf("There should be three corresponding files in the current folder.\n");
    printf("(1) <gene name>.fa, which is the DNA sequeence of the unspliced format of the gene.\n");
    printf("(2) <gene name>.bed, which contains the junctions informations.\n");
    printf("(3) <gene name>.sea, which contains the low coverage informations in wig format.\n");
    printf("The output will be the <gene name>_F.fasta and <gene name>_R.fasta\n");
    printf(", which is the enumearted isoforms in fasta format\n");
}

CJunctions getSplicingJunctions(string filename, bool forward, const char* pcaSeq = NULL)
{
    CJunctions js(filename, forward);
    getJunctionsFromBED(filename.c_str(), js, pcaSeq);
    js.sortJunctions(forward);
    printf("\rGot the junctions data");
    return(js);
}

CJunctions getLowCoverageRegions(string filename, bool forward)
{
    vector<pair <int, int> > oceans;
    const int COVERAGE_THRESHOLD = 0;
    getLowCoverageRegionFromWig(filename.c_str(), oceans, COVERAGE_THRESHOLD );
    printf("\rGot the low coverage region\n");
    CJunctions seas(filename, forward);
    seas.addJunctions(oceans, forward);
    seas.sortJunctions(forward);
    return(seas);
}

int enumerateIsoforms(int geneStart, int geneEnd, const char* contig, const char* geneName , bool forward)
{
    if (forward) {
        LOG_INFO("Info %d: Enumerate the forward isoforms\n", FINE_LOG);
    } else {
        LOG_INFO("Info %d: Enumerate the backward isoforms\n", FINE_LOG);
    }
    CJunctions js = getSplicingJunctions\
                    (string(geneName).append(string(".bed")), forward, contig);
    CJunctions seas = getLowCoverageRegions(string(geneName).append(string(".sea")), forward);
    // enumerate the isoform in the reverse strand
    CIsoformsEnumerator e(geneStart, geneEnd, forward, &js, &seas);
    int travelStart = forward ? geneStart : geneEnd;
    vector<int> possibleIsoformStarts = getPossibleIsoformStarts(travelStart, js, seas);
    for (vector<int>::iterator it = possibleIsoformStarts.begin();
            it != possibleIsoformStarts.end(); it++) {
        e.enumerateIsoforms(*it, MAX_EXONS_NO, e.geneEnd- e.geneStart);
    }
    cout << "There are " << e.isoforms.size() <<" isoforms."<< endl;
    // print only the coding seqeunces
    set<string> CDSs = getCDSs(e.isoforms, contig);
    string outFileName = string(geneName) + string("_CDSs.fa");
    int noOfCodingSeq = printCDSs(outFileName.c_str(), CDSs, string(geneName), forward);
    cout << "There are " << noOfCodingSeq << " CDS in " << outFileName << '.' << endl;
    // Alternatively, one can print the isoforms with long enough coding seqeunce

    vector<CIsoformInfo> selected_isoforms;
    set<string> mrnaSet;
    for (int i = 0; i < (int)e.isoforms.size(); i++) {
        char isoformId[FILENAME_MAX];
        char strand = e.forward ? 'F' : 'R';
        sprintf(isoformId, "I%c%d", strand, i);
        CIsoformInfo isoform = getIsoformInfo(geneName, isoformId, contig, e.isoforms.at(i));

        bool seqNotExist = (mrnaSet.find(isoform.seq) == mrnaSet.end());
        // filtering the isoform
        if (seqNotExist && (TEST || filterIsoform(isoform))) {
            selected_isoforms.push_back(isoform);
            mrnaSet.insert(isoform.seq);
        }
    }
    if ((int)selected_isoforms.size() > MAX_NO_ISOFORM_CANDIDATE ) {
        // Select a subset of isoform if there are too many candidates.
        selected_isoforms = selectByCDsScore(selected_isoforms, MAX_NO_ISOFORM_CANDIDATE);
    }
    LOG_INFO("Info %d: There are %d qualified isoforms in %s\n",\
             FINE_LOG, (int)selected_isoforms.size(), geneName);
    outFileName = string(geneName).append("_I.fa");
    printIsoforms(outFileName.c_str(), js.forward, selected_isoforms);
    return(0);
}

// print the isoform sequence as the fasta file
/*
int printIsoforms(ostream& out, vector<CIsoformInfo>& isoforms)
{
    char strand = this->forward ? 'F' : 'R';
    char filename[FILENAME_MAX];
    sprintf(filename, "%c_%s", strand, geneName.c_str());
    chExtName(filename,"isoform.exons");
	ofstream exonsFile(filename);
    int noQualifiedIsoform = 0;
    for (unsigned int i = 0; i < this->isoforms.size(); i++) {
        char isoformID[MAX_LINE];
        sprintf(isoformID, "I%c%d", strand, i);
		CIsoformInfo isoformInfo = getIsoformInfo(geneName, isoformID, contig, isoforms.at(i));
		if (TEST || filterIsoform(isoformInfo)) {
			reportIsoform(out, isoformInfo);
            exonsFile << "I" << strand << i << '\t' << isoformTag(isoforms.at(i)) << endl;
			noQualifiedIsoform ++;
        }
    }
	printf("\nThere are %d qualified isoforms in %s\n", noQualifiedIsoform, geneName.c_str());
	exonsFile.close();
	return(noQualifiedIsoform);

    return(0);
}*/

bool checkInputFileExist(const char* geneName)
{
    char filename[MAX_LINE];
    sprintf(filename, "%s.fa", geneName);
    if (!checkFileExist(filename)) return(false);
    sprintf(filename, "%s.bed", geneName);
    if (!checkFileExist(filename)) return(false);
    sprintf(filename, "%s.sea", geneName);
    if (!checkFileExist(filename)) return(false);
    return(true);
}

// This function set up some simple test suits
bool integralTest(void)
{
    cout << "testCIsoformsEnumerator: " << endl;
    bool pass = testCIsoformsEnumerator();
    cout << "testGetTheLongestORF: " << endl;
    pass = pass && testGetTheLongestORF();
    cout << "testPossibleStartSites: " << endl;
    pass = pass && testPossibleStartSites();
    if (pass) {
        cout << "OK" << endl;
    } else {
        cout << "FAIL" << endl;
    }
    return(pass);
}

int main(int argc, const char* argv[])
{
    if (argc != 2) {
        print_info();
        return(-1);
    }
    const char* geneName = argv[1];
    if (!checkInputFileExist(geneName)) {
        return(-2);
    }
    if (TEST) {
        integralTest();
    };

    char filename[MAX_LINE];
    sprintf(filename, "%s.fa", geneName);
    CchromosomeNTdata chr(filename);
    printf("\rGot the sequence data");
    int geneStart = 0;
    int geneEnd = chr.iChromosome_size - 1;
    enumerateIsoforms(geneStart, geneEnd, chr.caChromosome, geneName, true);
    enumerateIsoforms(geneStart, geneEnd, chr.caChromosome, geneName, false);
    return 0;
}
